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Abstract. A core-collapse supernova might produce large amplitude gravitational 
waves if, through the collapse process, the inner core can aquire enough rotational 
energy to become dynamically unstable. In this report I present the results of 3-D 
numerical simulations of core collapse supernovae. These simulations indicate that for 
some initial conditions the post -collapse inner core is indeed unstable. However, for 
the cases considered, the instability does not produce a large gravitational-wave signal. 
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^ : INTRODUCTION 

^ \ Core-collapse supernovae are potentially rich sources of gravitational waves. 

When a rotating stellar core exhausts its nuclear fuel the matter in the polar regions 
collapses more rapidly than the matter in the equatorial plane, since the latter must 
?h ! fight harder against centrifugal forces as it spirals inward. The changing oblateness 

of the core during this infall phase, and the subsequent changes due to core bounce, 
will generate gravitational waves (see, for example, Refs. [1,2]). After core bounce 
convective instabilities will cause the hot neutron star to boil [3] , and the resulting 
convective motions will give rise to a gravitational- wave signal [4]. To a lesser ex- 
tent, the same boiling process and release of gravitational waves can occur in the 
neutrino-heated material interior to the shock wave [4]. The supernova explosion 
itself may have a preferred direction, as evidenced by the large kick velocities of 
many neutron stars [5]. Such an asymmetry will generate a gravitational- wave 
signal with memory [6]. Perhaps the most interesting possibility, and the one dis- 
cussed in this paper, is that the stellar core will spin up as it collapses and produce 
a very rapidly rotating neutron star. The neutron star might be subject to dynam- 
ical instabilities that act to deform, or even fragment the star, and in the process 
produce large amplitude gravitational waves. 

In this paper I present the results of numerical simulations aimed at determining 
the types of initial conditions for a pre-collapse stellar core that lead to a dynam- 
ically unstable post-collapse inner core. The only previous investigations along 
these lines is found in the work of Rampp, Miiller, and Ruffert [7]. 



Analytical and numerical work on rapidly rotating fluid stars with Maclaurin- 
like rotational laws has shown that dynamical instabilities, in particular the m = 2 
bar-mode instability, will grow when the stability parameter T/|W| (the ratio of 
rotational kinetic energy to gravitational potential energy) exceeds about 0.27 [8]. 
For certain angular velocity profiles, a dynamical m = 1 instability can grow for 
T/|VF] as low as ~ 0.14 [9]. A back-of-the-envelope calculation frequently quoted 
in the literature suggests that the stability parameter will scale as T/|W| ~ 1/R 
during collapse, where R is the radius of the core. For a 10 solar mass star whose 
core collapses completely to neutron star size, this implies a factor ~ 100 increase 
in Based on this argument, one would expect that even the most slowly 

rotating cores will be dynamically unstable after collapse. 

The factor ~ 100 increase in T/|W| is overly optimistic for two reasons. First, 
centrifugal forces can halt the collapse at subnuclear density, producing a bloated 
inner core with below the threshold for dynamical instability. Second, the 

core actually "implodes" , rather than collapses, and some of the stellar core's matter 
and angular momentum remain outside the inner core. Since only a percentage of 
the core's angular momentum is drawn into the inner core, the increase in rotational 
kinetic energy is less than one would predict by assuming a more complete collapse. 
On the other hand, the simulations described here provide evidence that the entire 
core need not have T/|W| > 0.27 for the bar mode to grow on a relatively short 
timescale. It might be possible for the bar instability to grown as a secular process 
on a timescale of a few rotation periods, due to coupling between the inner and 
outer core regions. 

Two of the most rapidly rotating pre-collapse models considered here lead to 
dynamically unstable post-collapse inner cores. These models are unstable to the 
m = 2 bar mode, but the resulting bars have insufficient density and spatial extent 
to generate large gravitational-wave signals. In both cases, the contribution to 
the gravitational- wave signal from the bar deformation is smaller in amplitude, 
by a factor of ~ 2-5, than the purely axisymmetric signal generated by the core's 
changing oblateness during collapse and bounce. Overall, for sources in the Virgo 
cluster, these signals are too weak to be detected by laser interferometers operating 
in broadband mode. For galactic sources, these signals should be detectable. 

NUMERICAL CODE AND INITIAL MODELS 

The numerical code used for these investigations was written in collaboration 
with John Blondin. It uses Newtonian hydrodynamics and Newtonian gravity, ig- 
nores neutrino heating and cooling, and uses the relatively simple analytical equa- 
tion of state discussed by Zwerger and Miiller [2,7]. The gravitational- wave signal 
is computed in the quadrupole approximation. The hydro dynamical equations are 
solved with VH-1, written by Blondin, J. Hawley, G. Lindahl, and E. Lufkin. VH-1 
uses the piecewise parabolic method [10]. The Poisson equation for the gravita- 
tional potential is solved with multigrid techniques. 



Our code models the supernova in a minimal way, retaining just enough physics 
to capture the gravitational-wave signal in the leading-order quadrupole approx- 
imation. Nevertheless, the computation is a challenge due to the discrepancy in 
length scales involved. While the stellar core has a radius of 1000's of kilometers, 
the computational grid must have zone sizes of no more than ~ 1 km to support the 
steep density and velocity gradients in the inner core. A uniform 3-D Cartesian grid 
would require ~ 10 10 zones for these simulations, which is not feasible with current 
technology. One possible solution to this problem is to use a spherical-coordinate 
grid with non-uniform radial spacing, so the zones are smallest in the inner regions 
of the grid. The difficulty with this approach is that the zones become too narrow 
in the angular directions near the coordinate origin, and this drives the Courant- 
limited timestep to zero. The solution we have adopted is a nested grid scheme, in 
which the computational domain is covered by a sequence of Cartesian grids with 
increasing resolution and decreasing size. In this way only the inner-most region 
is covered with high resolution. This approach was also used by Rampp, Miiller, 
and Ruffert [7]. The simulations described here use 7 grids, each with 64 3 zones. 

In addition to the 3-D code just described, I also work with a 2-D code that 
assumes axisymmetry. The 2-D code uses a system of nested Cartesian grids in the 
r-z plane to achieve the necessary resolution of the inner core. The 2-D simulations 
reported here use 6 grids, each with 128 2 zones. 

In this paper I consider a sequence of initial models of the pre-collapse stellar 
core. The goal is to determine which of these models lead to dynamically unstable 
post-collapse inner cores. Each initial model is a rotating, equilibrium polytrope 
with T = 4/3 and central density p c = 10 10 g/cm 3 . The rotation laws for the 
sequence are given by 

n{zj) = fi e~ (ro/1500km)2 ; Q = 8 , 12, 16, 20, 24rad/s, (1) 

where Q is the angular velocity and w is the distance from the rotation axis. The 
Gaussian form for Q(w) was motivated in part by the results of Heger, Langer, and 
Woosley [11]. Their simulations of stellar evolution with rotation yield pre-collapse 
cores with relatively broad, Gaussian-like angular velocity profiles, similar in shape 
to the profile (1). Note, however, that even their most rapidly rotating model has 
an overall scale of Qq ~ 10, somewhat smaller than most of the models (1). 

In their previous work Rampp, Miiller, and Ruffert [7] used an initial data set 
with angular velocity profile Q(zu) = Oo/[1 + (cu/lOOkm) 2 ] and Q m 140rad/s. 
The resulting model is highly differentially rotating. The angular velocity drops 
from its central value of about 140 rad/s to less than 4 rad/s at a distance of R eq /2, 
where R eq ~ 1260 km is the equatorial radius. By contrast, the angular velocity 
of the most rapidly rotating model (1) drops from 24 rad/s to about 13 rad/s at 
R C J2, where R cq « 2500 km. 

The equation of state [2,7] contains a polytropic part whose stiffness depends 
on whether the matter is below or above nuclear density, and a thermal part that 
models the thermal pressure of matter heated by shock waves. Collapse of the 



initial models is induced by choosing the polytropic index in the sub-nuclear density 
regime to be T = 1.28, significantly below the value of 4/3 required for the model to 
remain in equilibrium. Random density perturbations at the 1% level were imposed 
at the beginning of the 3-D simulations. 



RESULTS 

The stability parameter is plotted as a function of time in Figure 1 for the 
sequence of models (1). These simulations were run with the 2-D code. Observe 
that T/|W^| increases by a factor of 2 or less as a result of core collapse. As 
shown in Figure 2, for the initial models with angular velocity Qq = 16, 20 and 
24 rad/s, the core experiences a centrifugal bounce and never reaches nuclear density 
Pnuc = 2.0 x 10 14 g/cm 3 . For the model with Qq = 12 rad/s, the inner core reaches 
nuclear density at core bounce then relaxes to a central value slightly below nuclear 
density. Only the most slowly rotating initial data, with Q = 8 rad/s, forms a stiff 
inner core with density greater than p nuc . These results show that centrifugal 
forces can severely inhibit core collapse, and prevent the stability parameter from 
experiencing the kind of growth suggested by the back-of-the-envelope arguments 
discussed in the introduction. 

The stability parameter for the two most rapidly rotating models exceeds 0.27 
after core bounce. I will use the labels f220 and £724 to denote the models with 
flo = 20 rad/s and Qq = 24 rad/s, respectively. For these models one might expect 
the post-bounce inner core, although centrifugally hung at sub-nuclear densities, 
to be dynamically unstable to growth of the m = 2 bar mode. The inner cores are 
not likely to be unstable to the m = 1 mode discussed in Reference [9], since their 
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FIGURE 1. T/\W\ vs. t for the sequence of models (1). 
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FIGURE 2. ln(p c ) vs. t for the sequence of models (1). 

density profiles are centrally peaked. Note, however, that prior to collapse both 
of these models exceed the nominal threshold T/|W| ~ 0.14 for growth of secular 
instabilities. Thus, these models are unrealistic — a realistic stellar core with such a 
high rate of rotation would lose its axisymmetry prior to collapse. With this caveat 
in mind, I have forged ahead and evolved these models with the 3-D code to check 
for dynamical instabilities in the post-collapse inner cores. These simulations show 
that model Q20 remains dynamically stable after collapse. Model f224 is unstable 
after collapse, and its inner core deforms into a bar shape. 

For the model considered by Rampp, Miiller, and Ruffert (RMR) [7], the pre- 
collapse core has an initial value for T/|W| of about 0.04, well below the threshold 
for growth of secular instabilities. Due to the high degree of differential rotation, the 
material in the outer layers rotates relatively slowly. As collapse proceeds, the lack 
of substantial centrifugal support allows the matter in the outer layers to strongly 
compress the inner core. The result is that the stability parameter increases by a 
larger fraction than that obtained with the data sets (1). The peak value of 
is about 0.35, although T"/|W| stays above 0.27 for less than 2 ms. After bounce, 
the inner core relaxes and T/|VF| quickly settles to a value of about 0.19. 

In their 3-D simulations, RMR imposed 10% random and 5% m = 3 density 
perturbations at a time of 2.5 ms before core bounce. The inner core showed 
m — 2, 3, and 4 asymmetries after bounce, but no significant enhancement in 
the gravitational- wave signal. Their simulations were halted at about 45 ms, ap- 
proximately 15 ms after core bounce. The following question naturally arises: are 
the post-bounce asymmetries seen in the RMR simulation merely a transient effect 
caused by the asymmetrical bounce? I have carried out a 3-D simulation with the 
RMR initial data using 1% random density perturbations imposed at the onset of 



collapse, to see if asymmetries will grow from a nearly axisymmetric bounce. I will 
refer to this simulation as model RMR. A second motivation for taking another 
look at the RMR initial data comes from the recent results in Reference [9], which 
suggest that a fluid body with toroidal density maximum (as occurs for model RMR 
both before and after collapse) can be dynamically unstable to an m = 1 mode for 
T/|VF] as low as ~ 0.14. The results of my simulation show that, as expected, 
the inner core is nearly axisymmetric immediately after core bounce. It remains 
axisymmetric until about 45 ms, which is the time at which RMR stopped their 
simulations. The dominant unstable mode that begins to grow at that time is the 
m = 2 bar mode, not the m — 1 mode. This occurs in spite of the fact that the 
stability parameter has a value of around 0.19. 

Before presenting the detailed results of the 3-D simulations, I need to establish 
some notation. The shape of the core can be described by expanding the matter 
density p in spherical harmonics: 

oo ra=l 

p(t, r, fl,0) = E E A im(t, r)Y em (0, 0) . (2) 

<=»m=-< 

The quadrupole formula relates the gravitational- wave amplitude to the second 
time derivative of the quadrupole moment of the mass distribution [12]. Inserting 
the expansion (2) into the quadrupole formula, we find 

C ^§h T + T = ] /- sin* Q(A 20 ) + y f (1 + cos* 0) n({A*)e»») 

— sine cose mA 21 )e^) , (3a) 
15 

C ^§h T J = -^coseS((i 22 )e^) - ^|sineS((i 21 >e^) , (3b) 

for the + and x components of the gravitational- wave amplitude in the transverse- 
traceless (TT) gauge. In these formulas R, 0, and $ specify the distance and 
angular direction from the source to the observation point and 3? and 3 denote real 
and imaginary parts. The angle brackets that appear in equations (3) are defined 
by (Ai m ) = / drr A Ap m , which is the spatial average of the second time derivative 
of A em , weighted with r 2 . 

From equations (3) we see that only the £ = 2 spherical harmonics contribute 
to the gravitational- wave signal in the quadrupole approximation. The coefficient 
A20 determines the oblateness of the mass distribution. Note that the space av- 
erage (A20) appears in h 7T , but not in h 7T . The coefficient A22 corresponds to 
a bar-shaped deformation in the equatorial plane. Growth of this coefficient im- 
plies growth of the usual Fourier m = 2 bar mode. The coefficient A 2 \ describes a 
bar-shaped deformation that is tilted out of the equatorial plane. This coefficient 
remains zero (apart from numerical noise) throughout the simulations, as one might 




expect from symmetry considerations. Note that the coefficient A u , which corre- 
sponds to the m — 1 mode of Reference [9], does not appear in the approximate 
formulas (3) for the gravitational- wave amplitude. 

Figures 3 and 4 show the ratios A20/P0 and |//Co for model f224, at various 
radii. Here, po is the average density at the given radius. For clarity of presentation, 
the results for radius r = 60 km are shown with heavy curves. From Figure 3 we see 
that, initially, the oblateness of the core increases rapidly as matter rushes inward 
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FIGURE 3. A 2 a/po vs. t for model 024. 
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FIGURE 4. \A 22 \/p vs. t for model 024. 



from all directions, but most rapidly along the polar regions. (A 2Q is negative for an 
oblate spheroid, positive for a prolate spheroid.) At ~ 28 ms the core experiences 
a "polar bounce" in which the matter in the polar regions is reflected off the inner 
core, reverses direction, and forms shock waves that propagate outward away from 
the equatorial plane. During the next ~ 10 ms, the oblateness of the inner core 
decreases as matter continues to rush in from the equator and out from the poles. 
At about 38 ms the core experiences an "equatorial bounce" in which the matter in 
the equatorial plane is reflected and forms an outwardly propagating shock wave. 
As the inner core relaxes, it spreads in the equatorial direction and again assumes 
a highly oblate shape. The oblateness remains fairly constant until near the end of 
the simulation, when the bar deformation becomes strong. Figure 4 shows growth 
of the bar mode, which begins around 100 ms. The spikes between 40 ms and 80 ms 
are caused by shock waves passing through the various radii. Since the shocks are 
not perfectly axisymmetric they produce relatively large but short-lived distortions 
that show up in the A 22 coefficient. The growth rate of the bar mode for model 
Q24 is d\n(\A 22 \/p )/dt « 46/s. 

A graph of the ratio |^422 1/ Po for model RMR shows that the bar mode begins 
to grow at about 45ms, at a rate of d\n(\A 22 \/ p Q )/dt ~ 180/s. | A 2 2 |/po reaches a 
peak value of 1.6 at ~ 70 ms. The strongest bar deformation occurs within a radius 
of ~ 100 km. 

The graphs in Figures 5 and 6 show the gravitational-wave signals for models 
Q24 and RMR, respectively. The solid curves show the + polarization amplitude 
(3a) as measured in the equatorial plane 6 = n/2 at a distance of R = 20 Mpc (the 
approximate distance to the Virgo cluster). The dashed curves show the + polar- 
ization amplitude (3a) as measured along the rotation axis = at R = 20 Mpc. 





Since the coefficient A 2 \ remains essentially zero, the x polarization amplitude (3b) 
at any angle is proportional to the + amplitude at = 0. Observe that for both 
models, the gravitational waves produced by the bar — the wiggles on the graphs at 
late times — are relatively small in amplitude compared to the gravitational waves 
produced by the core's (nearly axisymmetric) collapse and bounce. Fourier anal- 
ysis of the gravitational- wave signal for model f224 shows two peaks, one around 
5 Hz due to the collapse and bounce motion of the core, and the other around 
40 Hz due to the rotation of the bar-shaped inner core. The total energy radiated 
in gravitational waves for model f224, up to the end of the simulation, is only a 
few x 1CT 9 M c 2 . For model RMR, the gravitational- wave signal is spread across 
the frequency range 25-250 Hz. The total energy radiated for this model, up to the 
end of the simulation, is a few x 10 -7 M c 2 . 



DISCUSSION 

After collapse, models f224 and RMR are unstable to growth of the bar mode. 
On the other hand, model f220 is stable. The 3-D simulation of model fi20 was 
carried out to ~ 200 ms beyond core bounce, and the coefficient A 22 showed no 
signs of growth. At first sight these results might seem surprising: As shown in 
Figure 1, the stability parameter for models f224 and f220 after core bounce exceed 
the nominal threshold ~ 0.27 while the stability parameter for model RMR has 
a sustained, post-bounce value of less than 0.20. Of course, our understanding 
of the bar mode instability is based primarily on studies of isolated, equilibrium 
polytropes with Maclaurin-like rotation laws. There is little reason to believe that 
such a body would be a good approximation to a post-collapse stellar core and, 
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FIGURE 7. Density contours for model f224 at 199ms (left) and 229ms (right). The contour 
levels, in g/cm 3 , are 5.0 x 10 9 , 1.0 x 10 10 , 2.0 x 10 10 , 5.0 x 10 10 , 1.0 x 10 11 , and 5.0 x 10 11 . 

indeed, for the cases studied here, it is not. The value of the stability parameter for 
the post-collapse core is not a good diagnostic for the presence of the bar instability. 

Inspection of the data for models f224 and RMR shows that only the most dense 
regions of the post-collapse core participate in growth of the bar mode. For the 
model f224 in particular, the bar deformation is contained within the region with 
density p ^ 10 10 g/cm 3 . Figure 7 shows contour plots of the density in the equa- 
torial plane at 199 ms, when the bar has moderate strength, and at 229 ms, when 
the bar is near full strength. In the first plot, the contours below ~ 10 10 g/cm 3 
are nearly circular apart from some m = 4 noise caused by the Cartesian grid. 
By the time of the second plot, the bar-shaped region with p > 10 10 g/cm 3 has 
created "wakes" in the surrounding matter. These wakes form spiral arms that 
trail from the ends of the bar, and give rise to the spikes seen in the contour at 
5.0 x 10 9 g/cm 3 . Figure 8 shows the density as a function of radius in the equatorial 
plane for models f220 and fi24. For model f224 the density data is taken at 161 ms, 
near the beginning of the bar mode growth. Note that for both models, the density 
has a peak in the center then levels off to a value of about 10 10 g/cm 3 . Beyond 
~ 700 km, the density drops sharply. 

Motivated by the observations above, I will define the inner core for models 
f220 and f224 to be the region interior to ~ 10 10 g/cm 3 . Thus, we can view the 
post-collapse configuration as consisting of a dense inner core with p > 10 10 g/cm 3 
surrounded by relatively low density material. It is this inner core region that is 
unstable for model f224 and stable for model f220. The inner core for model RMR 
can be defined roughly by p > 10 10 g/cm 3 as well. 

Some insights into the behavior of the three models f220, Q24, and RMR can 
be gained by defining a stability parameter for the inner core, T ic /\Wi C \. The 
rotational kinetic energy of the inner core T ic is straightforward to compute. The 
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FIGURE 8. Density profile in the equatorial plane for model J124 at 161 ms and for model £120 
at 182 ms. 
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FIGURE 9. The stability parameter as a function of density cut-off for models 1124 at 161 ms, 
J120 at 182 ms, and RMR at 52.9 ms. 



gravitational potential energy Wi C must include the binding energy of the inner core 
material with itself, as well as the binding energy between the inner core material 
and the outer core material. The graph in Figure 9 shows the stability parameter 
of the region with density p > p cut , as a function of the cut-off value p cut , for the 
three models. The stability parameter for the inner core is obtained by setting 
p cu t ~ 10 10 g/cm 3 . The graph shows that for model Q20, T ic /\W ic \ ~ 0.15 while for 



models f224 and RMR, T ic /\W ic \ 0.19. Based on these results, we might expect 
that model RMR is less stable than model f220. Indeed, as the 3-D simulations 
show, RMR is unstable and Q20 is stable. Perhaps the more interesting question is 
this: Why are the inner cores for f224 and RMR unstable, given the fact that the 
values for their stability parameters are well below 0.27? The explanation might 
simply be that the post-collapse inner cores are not equilibrium polytropes with 
Maclaurin-like rotation laws, so the threshold for dynamical instability might be far 
different from 0.27. Another explanation might be that a post-collapse inner core 
is not isolated, and coupling to the outer core material can drive the bar instability. 
If this is correct, then growth of the bar mode is a secular process, as discussed 
for example by Schutz [13]. Note that for isolated, equilibrium polytropes with 
Maclaurin-like rotation laws, the threshold for growth of the secular instability is 
about 0.14. The inner core stability parameters for Q24 and RMR are well above 
this threshold. On the other hand, the stability parameter of the inner core for 
model f220 is very close to the threshold. 
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